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ABSTRACT 

The combined CfA2/SSRS redshift catalog has recently been improved and made public. We compute 
the redshift-space power spectrum of this catalog using the regions at —2.5° < (5i95o < 50° north 
{^^ < ai950 < 17'*) and south (20'* < aigso < 3'*, & < -13°) of the Galactic plane, where it is 98% 
complete down to Zwicky magnitude 15.5 and contains 13,681 galaxies. Our analysis uses Heavens- 
Taylor mode expansion, Karhunen-Loeve data compression and the Fisher matrix technique to compute 
quadratic band power estimates. This allows an exact calculation of window functions, including the 
integral constraint, in addition to the production of a power spectrum with uncorrelated error bars. Our 
results with this larger data set agree well with previous studies. We analyze 101/i~^Mpc and 130ft.~^Mpc 
volume-limited subsets in addition to the full magnitude-limited sample, and our results are well fit by, 
e.g., ACDM models with bias b —1.2, 1.4 and 1.4, respectively. We estimate the effect of extinction using 
the Schlegel, Finkbeiner & Davis dust map. Our results are exclusively for the redshift space galaxy 
power spectrum, so they can only be compared with theoretical predictions if appropriate corrections 
are made for biasing and redshift space distortions. 

Subject headings: large-scale structure of universe — galaxies: distances and redshifts — galaxies: 
statistics — methods: data analysis 



1. INTRODUCTION 

The three-dimensional maps of the Universe provided 
by galaxy redshift surveys provide valuable information 
about many fundamental cosmological parameters. This 
has motivated the creation of large data sets such as the 
Harvard-Smithsonian Center for Astrophysics (CfA), Las 
Campanas (LCRS) and IRAS 0.6 Jy (PSCz) redshift cata- 
logs, each well in excess of 10^ galaxies. Even more ambi- 
tious projects are currently under way, with the AAT two 
degree field survey (2dF) aiming for 250,000 galaxies and 
the Sloan Digital Sky Survey (SDSS) for 1 million. Com- 
bining such surveys with measurements from upcoming 
CMB experiments, e. 5., the MAP satellite, enables many 
cosmological parameters to be measured much more ac- 
curately than would be possible using CMB information 
alone (Eisenstein et al. 1999) — in principle. Achieving 
this in practice requires that the matter power spectrum 
P{k) can be accurately measured despite all the compli- 
cations that are inevitably present in real-world redshift 
surveys. 

This challenge has stimulated a large body of work over 
the last few years aimed at tackling such real- world issues, 
ranging from improved models of bias and extinction to 
the development of powerful new methods for measuring 
the power spectrum from surveys with arbitrary geometry 
and selection functions. 

The new CfA/SSRS UZC catalog has recently been com- 
pleted (Falco et al. 1999), made public and had its correla- 
tion function (Girardi et al. 2000) and topology (Schmalz- 
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ing & Diaferio 2000) measured. Since the power spectrum 
analyses of the original CfA2 and SSRS surveys were per- 
formed some time ago (Park et al. 1994; da Costa et al. 
1994; Marzke et al. 1994), it is therefore quite timely to 
apply some of these new methods to this extensive and 
further improved data set. This is the purpose of the 
present work. We limit our analysis to the galaxy red- 
shift space power spectrum on large scales, and therefore 
do not include the important complications of biasing, red- 
shift space distortions or nonlinearities. 

The rest of this paper is organized as follows. After de- 
scribing the UZC data set in Section |^, we present our 
analysis method in Section ^. The results are given in 
Section U, which also tests their robustness to various un- 
derlying assumptions. We summarize the conclusions and 
remaining challenges in Section |^. 

2. DATA 

The previous CfA redshift surveys (Huchra, Vogeley & 
Geller 1999; Huchra, Geller & Corwin 1995; Huchra et al. 
1990; Geller & Huchra 1989; Huchra et al. 1983; Davis 
et al. 1982) have all been based on the Zwicky catalog of 
galaxies with magnitude mzw < 15.5, involving a hetero- 
geneous sets of galaxy coordinates and redshifts. Falco 
et al. (1999) recently completed and made public the Up- 
dated Zwicky Catalog (UZC), which was improved over its 
predecessors in a number of ways: 

1. Over 5,000 redshifts have been re-measured or uni- 
formly re-reduced. 

2. The accuracy of the galaxy coordinates has been im- 
proved to better than two arcseconds using the dig- 
itized POSS plates. 



3. The redshift "blunder rate" 
substantially reduced. 



has been estimated and 



The sample described in Falco et al. (1999) consists of 
19285 galaxies, some without measured redshifts. The 
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Fig. 2.— The redshift dist/ih^itioMltSf'the 98% complete sample 
of 13,681 galaxies. The smooth curve is the selection function of de 
Lapparent et al 1989. 



sample has subsequently been further irnproved and ex- 
tended to 19415 galaxies, available onliner|. 

We use the latest (as of November 1999) version of this 
data set, which consists of 18,763 galaxies with redshifts. 
Their angular distribution is shown in Figure ^ Our anal- 
ysis is limited to the subset defined by 20'' < aigso < 4'' 
(hereafter "South") and 8'* < aigso < 17'* (hereafter 
"North"), -2.5° < ^1950 < 50°, where the sample is 98% 
complete down to Zwicky magnitude 15.5 (Falco et al. 
1999). We also apply a Galactic cut I&77I > 13° and dis- 
card the southern Galactic region with aigso > 3'* to re- 
duce extinction problems. This leaves a subset consisting 
of 13,681 galaxies in a region subtending about a quarter 
of the sky (10369 square degrees « 3.16 steradians). The 
redshift distribution of these galaxies is shown in Figure ^ 
Here and throughout this paper, we neglect the cosmic de- 
celeration correction and define the distance to a galaxy 
in redshift space as simply 



r= — = (3000/1^1 Mpc) z, 
H 



(1) 



since 2; <C 1, where z is the galaxy redshift in the CMB 
rest frame. Figure || also shows the radial selection func- 
tion that we use in our analysis, taken from de Lapparent 
et al. 1989. It assumes a Schecter luminosity function with 
parameters — —19.3 -I- 5 log h and a = —1.1. We trun- 
cate this magnitude limit radially so that 10~^ft,~^Mpc < 
r < 150/i~^Mpc, which leaves 13,184 galaxies. The lower 
limit removes 32 galaxies with negative redshifts. To fa- 
cilitate comparison with prior work, we also analyze two 
volume-limited subsamples of the data, where we keep only 
galaxies whose absolute magnitude is bright enough that 
they would have been visible above the mzw < 15.5 cut 
even if they were at the edge of the volume. Following 
Park et al. 1994, we take these to have depth 101 Mpc 
and 130 Mpc, which gives samples of size 2061 and 
909, respectively, with a selection function that should in 
principle be constant within the volume. 

3. METHOD 

® The latest version of UZC data set is available at 
cfa-www. harvard, edu/^ falco /UZC 



We analyze the data using the methods described in 
Tegmark et al. (1998, hereafter "T98"). The analysis con- 
sists of the following four steps: 

1. Heavens- Taylor pixclization 

2. Karhunen-Loeve compression 

3. Quadratic band-power estimation 

4. Fisher decorrelation 

We will now describe these steps in more detail. 

3.1. Step 1: Heavens- Taylor pixelization 

Our raw data consists of -/Vgai three-dimensional vec- 
tors Ta, a = 1, iVgai, giving the measured positions 
of each galaxy in redshift space. As discussed in T98, 
it is convenient to define the overdensity in "pixels" 
Xi, i = 1, ...,N^ by 



n(r) 
n{r) 



1 



n(r) 
n{r) 



il:,{Y)d?r ~ m, (2) 



for some set of functions tf^i and work with the iV^:- 
dimensional data vector x instead of the the 3 x iVgai num- 
bers . Since the mean density term 



has been subtracted out, we have 

(x)=0, 
(xxt) = C = N + S, 

where the shot noise covariance matrix given by 

_ r ^»(r)V.,(r) 3 
''~J n{v) 

and the signal covariance matrix is 



V;.(k)V^j(k)*P(fc) 



(2^)3 



(3) 



(4) 
(5) 



(6) 



(7) 



Here hats denote Fourier transforms and n is the three- 
dimensional selection function of the galaxy survey, i.e., 
n{Y)dV is the expected (not the observed) number of 
galaxies in a volume dV about r. 

We follow Heavens & Taylor (1995, hereafter "HT") in 
choosing our functions ipi to be spherical waves 



V'i(r) = Wi{r)Yirn(^)ji,{ki>„r), 



(8) 



where is a spherical harmonic and a spherical 
Bessel function. This choice of functions has the advantage 
(Fisher, Scharf, & Lahav 1994; HT) of forming an almost 
complete set, separating large and small scales fairly well 
and allowing efficient computation of m and C. The pix- 
elized data vector x is shown in Figure ^ for one of the 
volume-limited subsamples. 

Here and throughout, we use a single index i to refer to 
the triplet (imn) specifying an HT mode. The wavenum- 
bers kin are chosen as in HT so that the derivative of 



3 




w, 



to be of the form (Feldman, Kaiser & Peacock 1994) 



x" 




Fig. 3. — The triangles show the 1331 elements Xi of the data vec- 
tor X (the HT expansion coefficients) for the 130 h~^Mpc volume- 
limited sample. If there were no clustering in the survey, merely 
shot noise, they would have unit variance, and about 68% of them 
would be expected to lie within the blue/dark grey band. If our 
ACDM prior power spectrum were correct, then the standard devi- 
ation would be larger, as indicated by the shaded yellow/light grey 
band. The green/grey curve is the rms of the data points Xi, aver- 
aged in bands of width 25, and is seen to agree fairly well with the 
shaded band. 



jiiktnr) vanishes at a fixed radius r = i?niax and je{kinr) 
has n zeros before that. We choose the weight functions 



Wi{r) cx 



n(r) 



l + n(r)P(fcfa)' 



(9) 



normahzed such that the shot noise has unit variance, i.e., 
N,:,; = 1. Here P is our prior guess for the power spect rum, 



which we will discuss in more detail in Section 4.3.1. The 



results turn out to be rather insensitive to this weight func- 
tion choice — we will return to this issue in more detail in 



Section 4.3.3 



We avoid the complex issues described in Tadros et al. 
(1999) by using real- valued spherical harmonics, which are 
obtained from the standard spherical harmonics by replac- 
ing e™"^ by V2sinm(/), 1, v^cosmt/) for m < 0, m = 0, 
TO > respectively. HT show that apart from redshift- 
space distortions, the signal covariance matrix is given by 
S = ^P'S'*, where the diagonal matrix 

Pji' = tmn,l'm'n' = Sii'P{kin) (10) 

contains the effect of the power spectrum]^ and the matrix 



= /z..„(r),,(fc.„r)r,„(f),Hfc.'n-)r., 



which is ipi(j) expressed in (^TOn)-space and contains the 
relevant information about the survey geometry. When 
computing we use the elegant time-saving trick discov- 
ered by HT of expanding the angular part of the selection 
function in spherical harmonics and replacing the angular 
integral by a sum over Clebsch-Gordan coefhcients. 



^In fact, equation ( px[ ) 
the survey depth <^ -Rmax. 



in Section i.3.4 



is merely an approximation, valid when 
We will discuss this issue in more detail 
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Fig. 4. — The triangles show the 1166 elements yi of the 
compressed data vector y (the KL expansion coefficients) for the 
130 /i~^Mpc volume-limited sample. If there were no clustering in 
the survey, merely shot noise, they would have unit variance, and 
about 68% of them would be expected to lie within the blue/dark 
grey band. If our ACDM prior power spectrum were correct, then 
the standard deviation would be larger, as indicated by the shaded 
yellow/light grey band. The green/grey curve is the rms of the data 
points Xi, averaged in bands of width 25, and is seen to agree better 
with the yellow/light grey band than the blue/dark grey band. 



3.2. Step 2: Karhunen-Loeve compression 

As can be seen in Figure ^, most of the HT coefficients 
Xi have a fairly low signal-to- noise ratio, i.e., SufNu < 1 
and they are dominated by shot noise. Moreover, the HT 
modes form a complete set over the entire sky, so they are 
overcomplete for our case of the UZC catalog which cov- 
ers only a fraction of the sky. Both of these facts suggest 
that data compression may be possible, whereby almost 
all the cosmological information is retained in a smaller 
set of modes. Such compression would accelerate the ma- 
trix operations described in step 3, where the number of 
operations required grows as the number of modes cubed. 

We therefore subject the data vector x to Karhunen- 
Loeve compression. This method (Karhunen 1947) was 
first introduced into large-scale structure analysis by Vo- 
geley & Szalay (1996). It was recently applied to the 
Las Campanas redsift survey (Matsubara et al. 1999) and 
has been successfully applied to Cosmic Microwave Back- 
ground data as well, first by Bond (1995) and Bunn (1995). 
We define a new data vector 



B*x, 



(12) 



where b, the columns of the matrix B, are the Nx eigen- 
vectors of the generalized eigenvalue problem 

Sb = ANb, (13) 

sorted from highest to lowest eigenvalue A and normalized 
so that b^^Nb = I. This implies that 



which means that the transformed data values y have the 
desirable property of being uncorrelated. In the approxi- 
mation that the distribution function of x is a multivariate 
Gaussian, this also implies that they are statistically inde- 
pendent — then y is merely a vector of independent Gaus- 
sian random variables. Moreover, equation (|l^) shows that 
the eigenvalues Xi can be interpreted as a signal-to-noise 
ratio S/N . Since the matrix B is invertible, the final data 
set y clearly retains all the information that was present in 
X. In summary, the KL transformation partitions the in- 
formation content of the original data set x into chunks 
that are 

1. mutually exclusive (independent), 

2. collectively exhaustive (jointly retaining all the in- 
formation), and 

3. sorted from best to worst in terms of their informa- 
tion content. 

Figure ^ shows that most of our KL coefficients yi have a 
signal-to- noise ratio A <C 1, so that the bulk of the cosmo- 
logical information is retained in the first Ny coefficients, 

Ny < Nx. 

In our case, numerical considerations force an additional 
compression step. The above-mentioned overcompleteness 
means that certain linear combinations of modes vanish 
almost completely within the angular mask of the survey 
(shown in Figure |l]). This makes many of of the eigen- 
values of S and N essentially zero and poses a numerical 
problem when solving equation (p^, since the standard 
reduction to an ordinary eigenvalue problem by Cholesky 
decomposing N or S fails when both are singular. Adding 
a small number to the diagonal elements of N also fails to 
solve the problem: since these redundant modes are tiny, 
avoiding to vanish completely mostly because of rounding 
errors, they have a minute contribution from both signal 
and noise and therefore will not automatically get weeded 
out by a cut on the signal-to-noise A. We therefore begin 
our compression by expanding the data in the eigenvectors 
of S and get rid of these redundant junk modes by throw- 
ing away all eigenmodes with eigenvalues below 10~®. We 
then subject the new compressed data set to KL compres- 
sion. 

In conclusion, this step takes the vector x and its co- 
variance matrix C from Figure |^ and compresses it into 
the smaller vector B*x and its covariancc matrix B*CB, 
illustrated in Figure ^. 

3.3. Step 3: Quadratic band-power estimation 

In this step, we perform a much more radical data com- 
pression by taking certain quadratic combinations of the 
data vector that can easily be converted into power spec- 
trum measurements. 

We parametrize the power spectrum as a piecewise con- 
stant function with Np "steps" of height pi , which we term 
the band powers and group into an TVp-dimensional vec- 
tor p. Thus P{k) = Pi for ki-i < k < ki, where 



(14) 



= ko <ki < ... < kNp. (15) 

We use Np — 52 bands, most of them with bandwidth 
A/c = 0.01//i~^ Mpc, which should provide fine enough k- 
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Fig. 5. — The triangles show the 52 elements qi of the raw 
quadratic band-power estimators q for the 130 /i~^Mpc volume- 
limited sample. If the fiducial ACDM power spectrum were correct, 
then about 68% of them would be expected to lie within the shaded 
yellow/light grey band, centered on the solid curve. 



resolution to resolve any baryonic wiggles and other spec- 
tral features that may be present in the power spectrum. 
This means that we can write 



4=0 



(16) 



where the derivative matrix C, = dC/dpi is the contribu- 
tion from the i*'* band and is computed by simply limiting 
the implicit sums in the matrix multiplication S = \I>PvI>* 
to those lines and columns of the power spectrum ma- 
trix P where lies in the i*'* band; ki-i < kin < ^i- 
For notational convenience, we included the noise term 
in equation ( p^ ) by defining C,o = N, corresponding to 
an extra dummy parameter po = 1 giving the shot noise 
normalization. 

Our quadratic band power estimates are defined by 



(17) 



i — 1, ...,Np. These numbers are shown in Figure ^ and 
we will group them together in an A'p-dimensional vec- 
tor q. Note that whereas x and y (and therefore C and 
B) were dimensionless, p has units of power, i.e., volume. 
Equation (^^ therefore shows that q has units of inverse 
power, i.e., inverse volume. It is not immediately obvious 
that the vector q is a useful quantity. It is certainly not 
the final result (the power spectrum estimates) that we 
want, since it does not even have the right units. Rather, 
like y, it is a useful intermediate step. In the approxi- 
mation that the pixelized data has a Gaussian probability 
distribution (a good approximation in our case because of 
the central limit theorem, since iVgai is large) q has been 



shown to retain all the information about the power spec- 
trum from the original data set (Tegmark 1997, hereafter 
"T97"). The numbers qi have the additional advantage 
(as compared with, e.g., maximum-likelihood estimators) 
that their properties are easy to compute: their mean and 
covariance are given by 



(q) 
(q>(q>* 



Fp, 



(18) 
(19) 



where F is the Fisher information matrix (Tegmark et al. 
1997) 

(20) 



F,, = -tr [C-iC^,C-iC,,] 



Quadratic estimators were first derived for galaxy survey 
applications (Hamilton 1997ab). They were accelerated 
and first applied to CMB analysis (T97; Bond, Jaffe & 
Knox 2000). 

In conclusion, this step takes the vector y and its co- 
variance matrix B*CB from Figure ^ and compresses it 
into the smaller vector q and its covariance matrix F, il- 
lustrated in Figure ||. Although equation (|l^) shows that 
we can obtain unbiased estimates of the true powers p by 
computing F~^q, there are even better options, as will be 
described in the next subsection. 

3.4. Step 4-' Fisher decorrelation 

Let us first eliminate the shot-noise dummy parameter 
Pq, since we know its value. We define f to be the 0*'* 
column of the Fisher matrix defined above {fi = F^o) and 
restrict the indices i and j to run from 1 to Np from now 
on, so f , q and p are A'p-dimensional vectors and F is an 
Np X Np matrix. Since po = 1, equation ( p^ ) then becomes 
(q) = Fp + f . 

We now define a vector of shot noise corrected band 
power estimates 

p = M(q-f), (21) 

where M is some matrix whose rows are normalized so 
that the rows of MF sum to unity. Using equations (18) 
and (n9h, this gives 



(PP* 



(P> 
(P>(P>* 



Wp, 
MFM*, 



(22) 
(23) 



where W = MF. We will refer to the rows of W as win- 
dow functions, since they sum to unity and equation (22) 
shows that pi probes a weighted average of the true band 
powers pj, the i*'' row of W giving the weights. 

3.4.1. The minimum-variance choice 

What is the best choice of the matrix M? A simple and 
natural choice is 



(24) 



I.e., M diagonal (T97; Bond, Jaffe & Knox 2000). Figure| 
plots the rows of the Fisher matrix, which have the same 
shape as the window functions for this choice of M. As is 
seen in this figure, all elements of F are positive. This is 
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Fig. 6. — The rows of the Fisher matrix are shown for the 
130 /i~^Mpc volume-limited sample. The i*'' row typically peaks 
at the i**" band, the scale k that the band power estimator qi was 
designed to probe. Large amplitudes signify large information — if 
all curves had exactly the same shape, then the area under the i*'' 
curve would be (Api)~^, the inverse variance of qi when normal- 
ized as a band power estimator. The turnover in the envelope at 
k ~ 0.3/h-^Mpc reflects our omission of HT modes probing smaller 
scales. 



guaranteed by the fact that the matrices in equation ( |20[ ) 
are ah positive definite^. 



3.4.2. The unbiased choice 
Another interesting choice is (T97) 

M ^ F-\ 



(25) 



which gives W = I. In other words, all window functions 
are Kronecker delta functions, and p gives completely unbi- 
ased estimates of the band powers, with (p) = pi regardless 
of what values the other band powers take. A drawback 
of this choice is that the new covariance matrix of equa- 
tion ( p^ ) becomes F~^, which usually gives substantially 
larger error bars (Ap.^ = M|/^ — [(F^^)^]^/^) than the 
first method, anti-correlated between neighboring bands. 

3.4.3. The uncorrelated choice 

The two above-mentioned choices for M both tend to 
produce correlations between the band power error bars. 
The minimum- variance choice generally gives positive cor- 
relations, since the Fisher matrix cannot have negative el- 
ements, whereas the unbiased choice tends to give anticor- 
relation between neighboring bands. The choice (Tegmark 
& Hamilton 1998; Hamihon & Tegmark 2000) 



E 



-1/2 



/2 



(26) 



* This can be seen as follows. F; 



UrfAjA,!, where A; 



C" i/2q^^q— 1/2 'pjjg trace of a such a product of two matrices is 
positive if both matrices are positive definite — this is obvious in the 
basis where one of them is diagonal, since neither can have negative 
diagonal elements. C^i is positive definite since it is a covariance 
matrix (corresponding to a power spectrum equal to unity in the i^^ 
fc-band and vanishing elsewhere), so the matrices are all positive 
definite. 



Fig. 7. — The window functions are shown for the 130 /i~^Mpc 
volume-limited sample using the decorrelation method. The i*^ row 
of W typically peaks at the i^^ band, the scale k that the band power 
estimator pi was designed to probe. To facilitate comparison with 
Figure tl we have rescaled the window functions to sum to (Ap;)"'^ 
rather than unity here. 



i.e., M = with the rows renormalized, has the at- 

tractive property of making the errors uncorrelated, with 
the covariance matrix of equation (23) diagonal. The cor- 
responding window functions W are plotted in Figure ^ 
and are seen to be quite well-behaved, even narrower than 
those in Figure ^ while remaining positive. This choice is 
a compromise between the two first ones: it narrows the 
minimum variance window functions at the cost of only a 
small noise increase, with uncorrelated noise as an extra 
bonus. Note that all three choices of M retain all the cos- 
mological information, since the vector y can always be re- 
covered by multiplying p by M~^. The minimum- variance 
band power estimators are essentially a smoothed version 
of the uncorrelated ones, and their lower variance was paid 
for by correlations which reduced the effective number of 
independent measurements. 



3.5. Integral constraint correction 

Since the main focus of this paper is the power spec- 
trum on the largest scales, it is crucial that we deal with 
the complication known as the integral constraint (Pea- 
cock & Nicholson 1991). If we knew the selection function 
n(r) a priori^ before counting the galaxies in our survey, 
we would be able to measure the power on the scale of the 
survey. Our power spectrum estimate would essentially be 
the square of the ratio of the observed and expected num- 
ber of galaxies in our sample. Of course, we do not know 
n a priori^ so we use the galaxies themselves to normalize 
the selection function. Thus the measured density fluc- 
tuation automatically vanishes on the scale of the survey, 
and and nai've application of the power spectrum estima- 
tion method we have described will falsely indicate that 
P{k) ^ as fc — + 0, regardless of the behavior of the true 
power spectrum on large scales. 

Fortunately, this problem has a simple remedy once the 
data has been pixelized. As described in T98, the problem 
is that we do not know the true amplitude of the mean 
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Fig. 8. — Decorrelated band-power measurements are shown for 
the 101/i~^Mpc volume-hmited sample. Although the window func- 
tions overlap slightly, the errors are all uncorrelated. The error bars 
include the effects of both shot noise and sample variance corre- 
sponding to the fiducial model. The horizontal bars show the rms 
width of the window function corresponding to each measurement. 
If the fiducial model (heavy curve) is correct, the measurements 
should on average equal this curve convolved with the window func- 
tions (dotted curve). The thin red line shows the shot noise contri- 
bution that has been subtracted out. 



density mode m defined in equation (^), and consequently 
may not have subtracted it out correctly in equation (||). 
We can therefore immunize our data to this problem by 
making it orthogonal to m. Defining a projection matrix 



Fig. 9. — Same as previous figure, but for the 130/i ^ Mpc volume- 
limited sample. 



n 



mf 



(27) 



the recipe is simply to replace x by IIx and the matrices 
C,i by TlC,iTl (T98). This trick of choosing modes that 
are orthogonal to the mean density was first suggested by 
Fisher et al. (1993). We find that this correction increases 
our error bars slightly on the largest scales, mainly in the 
leftmost band. 

4. RESULTS 

4.1. Basic results 

Our basic results are shown in Figures |, I and |0|. 
Here the band powe r esti mates are computed using the 
M-choice of Section |3.4.3| , so the error bars are uncorre- 
lated. The horizontal location of each data point and the 
width of the horizontal bars are the mean and the stan- 
dard deviation of the corresponding window function, re- 
spectively. To avoid excessive clutter, we have averaged 
neighboring band-power measurements (and their corre- 
sponding window functions) on the smallest scales using a 
simple inverse-variance weighting since the error bars are 
uncorrelated. 

The heavy curve shows the prior power spectrum vec- 
tor p that was used in the calculation. It is a flat ACDM 




0.01 
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Fig. 10. — Same as previous figure, but for the magnitude limited 
sample. 



"concordance model" (Wang et al. 1999) with Q\ = 0.7, 
^cdm = 0.25, fif, = 0.05, h = 0.65, n = 1 with tr = 1 for the 
matter fluctuations, rescaled by a different bias factor for 
each of the three galaxy samples. The linear power spec- 
trum for this model was computed using the fit of Eisen- 
stein & Hu (1999), then corrected for nonlinear effects with 
the formalism of Jain et al. (1995) based on HKLM scaling 
(Hamilton et al. 1991; Peacock & Dodds 1996), with the 
local power spectrum slope given by the "baryon wiggle- 
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free" version of the spectrum. Our non-linear as — 1 nor- 
malization corresponds to a linear erg « 0.93. We first 
estimated a power spectrum assuming 6=1, then iter- 
ated the calculation once with new bias factors provid- 
ing a better normalization to the actual measurements. 
For the 101 /i~^Mpc, 130/i~^Mpc and magnitude-limited 
samples, these three bias factors b (reflected by the heavy 
curves shown in Figures H, ^ and 10) are 1.2, 1.4 and 1.4, 
respectively. 

This is not a complete treatment of non-linearity (see, 
e.g., Meiksin & White 1999; Scoccimarro et al. 1999; 
Hamilton 2000; Benson et al. 1999). The effects of non- 
linearity will not bias our power spectrum estimates, but 
the error bars are likely to be underestimated on nonlinear 
scales k > 0.3//i~iMpc. 

Since each data point is the power spectrum convolved 
with a window function, we would not expect the data 
points to fall exactly on the true power spectrum even on 
average. Rather, if the prior were correct, they would on 
average fall on the dotted curve Wp, which shows the prior 
averaged with the window functions of each band. In all 
cases, the results are seen to be consistent with the priors 
used, except perhaps on the very smallest scales and for the 
leftmost band in the ma gnitude-limited case, to which we 
will return in Section 4.S below. We will discuss robustness 



towards the choice of prior below in Section 4.3.1. 

How reliable are these results? In the remainder of this 
section, we present a series of tests, both of our software 
and algorithms and of potential systematic errors. 

4.2. Validation of method and software 

Since our analysis consists of a number of somewhat 
complicated steps, it is important to test both the soft- 
ware and the underlying methods. We do this by gener- 
ating A'montc — 200 Monte Carlo simulations of the UZC 
catalog with a known power spectrum, processing them 
through our analysis pipeline and checking whether they 
give the correct answer on average and with a scatter cor- 
responding to the predicted error bars. We found this 
end-to-end testing to be quite useful in all phases of this 
project — indeed, things worked on neither the first nor 
the second attempt... 

4.2.1. The mock survey generator 

Standard N-body simulations would not suffice for our 
precision test, because of a slight catch-22 situation: the 
true non-linear power spectrum of which an N-body simu- 
lation is a realization (with shot noise added) is not known 
analytically, and is usually estimated by measuring it from 
the simulation — but this is precisely the step that we wish 
to test. We therefore resort to a simpler approach, where 
we generate realizations that are still in the linear regime, 
just as is routinely done when preparing initial conditions 
for N-body simulations. We do this in the following steps: 

1. Generate a Gaussian random field with the pre- 
scribed power spectrum on a cubic grid (by gener- 
ating uncorrelated Gaussian random variables with 
variance P{k) at each grid point in Fourier space, 
then performing an FFT). 

2. Generate a random galaxy position inside the an- 
gular mask with radial probability distribution pre- 
scribed by the selection function. 




Fig. 11. — The triangles show the elements Xi of the data vector 
X (the HT expansion coefficients) averaged over 200 Monte-Carlo 
simulations of the 101 /i~^Mpc volume-limited sample. If the algo- 
rithms and software are correct, then their mean should be zero and 
about 68% of them should lie within the shaded yellow/grey region 
giving their standard deviation. The green/grey curve is the rms of 
the data points Xi, averaged in bands of width 25, and is seen to 
agree well with a smoothed version of the shaded region. 



3. Evaluate the density field S at using trilinear in- 
terpolation between the nearest points on the FFT 
grid. 

4. Add this galaxy to the catalog with a probability 
equal to [1 -I- (5(ri)]/2, otherwise discard it. 

5. Go back to step 2 and repeat until the desired num- 
ber of galaxies has been generated. 

To avoid the quantity 1 -I- S{ri) going negative in step 4, 
which would spoil our procedure, we normalize our fidu- 
cial power spectrum to have much less power than ob- 
served in the actual Universe. We choose our test power 
spectrum to be a simple Gaussian P{k) oc e~'^^^'^ with 
R — 32/i^^Mpc, normalized so that the rms fluctuations 
((52)1/2 ^ Q 2. This ensures that \5\ > 1 (breaking step 4) 
occurs only a negligible fraction of the time (about once 
in 1.7 million). 

4.2.2. Testing the Heavens- Taylor pixelization 

Figure |ll| shows the result of processing the Monte Carlo 
simulations through the first step of the analysis pipeline, 
i.e., computing the corresponding Heavens- Taylor expan- 
sion coefficients Xi. This is a very sensitive test of the 
mean correction given by equation (^, which can be a 
couple of orders of magnitude larger than the scatter in 
Figure |ll| for some modes. A number of problems with 
the radial selection function integration and the spherical 
harmonic expansion of the angular mask in our code were 
discovered in this way. After fixing these problems, the co- 
efl[icients Xi became consistent with having zero mean as 
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Fig. 12. — The triangles show the elements qi of the raw quadratic 
band-power estimators q, averaged over 200 Monte-Carlo simula- 
tions of the 101 /i~^Mpc volume-limited sample. If the algorithm 
and the software is correct, then about 68% of them would be ex- 
pected to lie within the shaded yellow/grey band, centered on the 
solid curve. 

seen in the figure. The figure also shows that the scatter 
in the modes is consistent with the predicted standard de- 
viation (Ji = {Cii/NjnontcY^^ (shaded region), with most 
of the the fluctuations being locahzed to modes probing 
large scales (with i, m and n being small). 

Processing the Monte Carlo simulations through the sec- 
ond step of the analysis pipeline showed that the corre- 
sponding KL eigenmodes yi passed the same test. 

4.2.3. Testing the quadratic band-power estimation 

Figure [ij shows the result of processing the Monte 
Carlo simulations through the third step of the analysis 
pipeline, i.e., computing the raw unnormalizcd quadratic 
band-power estimates qi. Since information from large 
numbers of modes contributes to each qi, the scatter 
(Ji = (Fii /Almonte) is seen to be small. Therefore, even 
quite subtle bugs and inaccuracies can be (and were!) dis- 
covered and remedied as a result of this test. 

4.2.4. Testing the Fisher decorrelation 

Figures ^ and |l^ show the result of processing the 
Monte Carlo simulations through the fourth and final step 
of the analysis pipeline, i.e., computing the decorrelated 
and normalized band-power estimates pi . The mean recov- 
ered power spectrum is seen to be in excellent agreement 
with the Gaussian prior used in the simulations (Figure 
convolved with the window functions, and the observed 
scatter is seen to be consistent with the predicted error 
bars (Figure M). These two fi gures therefore constitute 
an end-to-endtest of our data analysis pipeline, since er- 
rors in any of the many intermediate steps would have 
shown up here at some level. 



1.5x10" - ' ' ' ' 




0.1 O.a 0.3 0,4 0,5 

k [1/h-' Mpc] 



Fig. 13. — The triangles show the decorrelated band-power 
estimates pi, averaged over 200 Monte-Carlo simulations of the 
magnitude-limited sample. If the algorithm and the software is 
correct, then about 68% of them should lie within the shaded yel- 
low/grey band, centered on the Gaussian fiducial power spectrum 
(solid curve). 




k [1/h ' Mpc] 



Fig. 14. — Same as the previous figure, but for the error bars 
Api. The triangles show the observed scatter in the 200 simulations. 
If the algorithm and the software is correct, then about 68% of 
them should lie within the shaded yellow/grey band predicted by 
the Fisher matrix formalism, centered on the solid curve. 

4.3. Robustness to method details 

Our analysis pipeline has a number of "knobs" that can 
be set in more than one way. This section discusses the 
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Fig. 15. — Effect of prior. The triangles show the decorrelated 
band-power estimates pi, averaged over 200 Monte-Carlo simula- 
tions of the magnitude-limited sample with a Gaussian power spec- 
trum (red line). The power spectrum estimation was carried out 
assuming a totally different fiducial model, the ACDM model to the 
right. If the data was consistent with the prior, about 68% of the 
measurements should fall within the narrow yellow/grey region. In- 
stead, the method is seen to faithfully reproduce the actual input 
power spectrum, with no evidence of a bias towards the assumed 
prior. 



sensitivity to such settings. 

4.3.1. Effect of changing the prior 

The analysis method employed assumes a "prior" power 
spectrum via equation (p^, both to compute band power 
error bars and to find the galaxy pair weighting that min- 
imizes them. As mentioned, an iterative approach was 
adopted starting with a simple ACDM model with cts = 1, 
then rescaling it to better fit the resulting measurements 
and recomputing the measurements a second time. To 
what extent does this choice of prior afi'ect the results? 
On purely theoretical grounds (e.g., Tegmark, Taylor & 
Heavens 1997), one expects a grossly incorrect prior to give 
unbiased results but with unnecessarily large variance. If 
the prior is too high, the sample-variance contribution to 
error bars will be overestimated and vice versa. This hy- 
pothesis has been extensively tested and confirmed in the 
context of power spectrum measurements from the Cosmic 
Microwave Background {e.g., Bunn 1995). An analogous 
test is shown in Figure Eq, showing that the correct result 
is recovered even when our 200 simulations are analyzed 
with a grossly incorrect prior. 

Generally, the pair weighting strives to minimize the 
joint contribution from sample variance and shot noise to 
the scatter in the measurements. This scatter will there- 
fore be unnecessarily large both if the prior is too low (so 
that sample variance not taken seriously enough) and if it 
is too high (so that excessive paranoia about sample vari- 
ance gives a pair weighting producing unnecessarily large 



105 




0.01 0,1 
k [1/h ' .Mpc] 

Fig. 16. — Same as Figure hcj (magnitude-limited sample), but 
for the method giving maximally narrow window functions. 
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Fig. 17. — Same as Figure hcj (magnitude-limited sample), but 
for the method giving smaller, oandwidth independent error bars. 



shot noise). The resulting scatter therefore increases only 
to second order when the prior is slightly off. Since the 
iterated priors used in our analysis of the real data agree 
so well with the actual measurements, slight remaining de- 
viations of the prior from the truth are therefore likely to 
have a negligible impact. 

4.3.2. Effect of changing the decorrelation method 
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Our main results for the power spectra (Figures p|-^10|) 
were com put ed using the decorrelation method given by 
equation (|2(^). To assess the sensitivity to this choice, we 
repeated the analysis for the other methods. The results 
were generally consistent with the same fiducial models, 
but as expected, the nature of the scatter was found to 
be strongly method-dependent. This is illustrated in Fig- 
ures 0and|l^ for the magnitude-limited case. 

Figure ^ used the method given by equation (^) , giv- 
ing maximally narrow window functions. Although they 
are plotted as having zero width, a calculation with nar- 
rower bands would show them to have a width of order that 
of the bands used here, i.e., of order the horizontal sepa- 
ration between neighboring points on the plot. However, 
Figure ^ also shows that there are no free lunches: the 
errors are bigger than in Figure |l^, and the error covari- 
ance matrix shows that they are strongly anti-correlated 
with their nearest neighbors. This figure is essentially a 
window-deconvolved version of Figure and smoothing 
it would recover that figure. 

Figure |l^ used the method given by equation (|2j) , and 
deviates from Figure |l^ in the opposite way from the pre- 
vious example. It is essentially a smoothed version of 
Figure giving nice small error bars, slightly broader 
window functions and positively correlated errors between 
neighboring points. The broader window functions are 
seen to be particularly annoying on the largest scales, 
where they shift the effective wavenumber probed far to 
the right. 

4.3.3. Effect of changing the galaxy weighting 

When expanding our galaxies in HT modes, we applied 
the radial FKP weighting given by equation (||). How does 
this particular choice affect the results? To address this 
issue, we repeated the analysis with equation (^) replaced 
by w{r) oc n{r), i.e., the simple P{k) = limit of equal 
weight per galaxy. This resulted in almost no perceptible 
loss of information, typically increasing the band-power 
error bars on large scales by less than a percent. The 
resulting power spectrum measurements were essentially 
unchanged. 

In the limit where infinitely many HT modes would be 
used, any functions whatsoever can be created by taking 
linear combinations of them, since they form a complete 
basis over the survey volume. This would make the choice 
of the radial weighting function w(r) completely irrelevant, 
since the subsequent KL compression step would end up 
recovering the true KL eigenfunctions regardless. The rea- 
son that the radial weighting makes any difference at all 
in our case is therefore that we have used only a limited 
number of modes to start with, making it important that 
they do not grossly over- or underweight sparse distant 
galaxies relative to nearby ones. 

4.3.4. Effect of other method details 

The analysis above was performed using 11"^ = 1331 
HT modes, with an angular cutoff at ^max — 10 giving 
(1 -|- i'max)^ = 121 angular modes Yfm and a radial cutoff 
"-max = 10 giving 11 radial modes jt{kinr) per Yn^. To 
explore the sensitivity to these choices, we repeated the en- 
tire analysis with £„iax = ?^inax = 3, 5 and 7 giving 64, 216 
and 512 HT modes, respectively. As expected, we found 
that the band powers on the very largest scales converged 
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Fig. 18. — Our 52 k-bands are shown (top and bottom) together 
with the discrete spectra of wavenumbers fcf „ used in the HT for- 
malism. Of the two spectra, the upper one corresponds to the values 
used in equation (H) for expanding the data, whereas the lower 
one corresponds to the values used in equation ( |lO| ) for comput- 
ing the covariance matrix C and its derivatives CTT, both for the 
magnitude limited analysis. 

quite rapidly as more modes were added, and that the 
new information made a difference mainly on the smaller 
scales where the new modes were sensitive. These numer- 
ical experiments suggested that our 1331 mode analysis 
retained a large fraction of the cosmological information 
down to fc ~ 0.3, and that a more ambitious analysis with 
more modes would give substantially smaller error bars at 
smaller scales. The information matrix plotted in Figure]^ 
reenforces this conclusion. 

To determine the best tradeoff between angular and ra- 
dial modes, we performed a number of tests with fmax 7^ 
n-max, keeping the total number of modes (1 -I- ^max)^(l + 
"■max) roughly constant. These tests indicated that the 
largest Fisher information (the smallest error bars) on the 
power spectrum on the largest scales was obtained for 
the symmetric case ^max ~ "max- This is why we used 
^max = "max in our main analysis. 

We used 52 fc-bands as illustrated in Figure [l^. To com- 
pute the covariance matrix C and its 52 derivatives C,i 
exactly using equation (|o|) , the contributions from a very 
large set of fc^m-values would be required. In practice, we 
truncated the fc^m-values used for these internal computa- 
tions (the lower subset in Figure |lj) at ^ = €max = 21, 
since experimentation with different cutoff values ^max 
showed that all results for fc < 0.3 had converged by then. 
An alternative approach to truncation is to make a sharp 
cut at a fixed fc-value (HT). 

As mentioned above, equation ( |l0|) is only approximate. 
The discrete spectra shown in Figure |l^ would apply only 
if the radial Bessel functions were allowed to extend to in- 
finite radius. Since they are truncated outside the survey 
volume, the sharp lines of Figure |l^ get slightly smeared 
out: they essentially get convolved with the Fourier trans- 
form of the survey volume, which gives them a charac- 
teristic width Afc of order the inverse size of the survey 
volume. 

In the eigenmode expansion, the cutoff was placed at 
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an S-eigenvalue of 10 ®, simply to avoid numerical sin- 
gularities. This reduced the original 1331 modes to 1182, 
1166 and 1206 for the 101/i"^Mpc, 130/i"^Mpc and mag- 
nitude limited samples, respectively. The results remained 
completely unchanged when this 10~^ cutoff was increased 
by two orders of magnitude. Since 1331 models are very 
manageable numerically, our principal motivation for per- 
forming the subsequent S/N-analysis was to check for sys- 
tematic errors. We therefore discarded no further modes 
in this step. 

4.4. Robustness to extinction 

Mis-estimates of n constitute a potential source of sys- 
tematic errors. Although our method for dealing with the 
integral constraint immunized against errors in the overall 
normalization of n, errors in its shape would still add spu- 
rious power to our estimate of P{k). Let us first consider 
errors in the angular part of n caused by Galactic dust, re- 
turning to errors in the radial part in the next subsection. 

The angular modulations caused by dust extinction 
tend to have a power spectrum rising sharply towards the 
largest scales (Vogeley 1998), and is therefore of particular 
concern for the interpretation of our leftmost bandpower 
estimates. Although the UZC subset we are using is 98% 
complete relative to the underlying Zwicky sample, extinc- 
tion can of course cause completeness modulations in the 
latter. To estimate the severity of this effect, we used the 
recent extinction map produced by Schlegel, Finkbeiner 
& Davis (1999), using their B-band conversion factor of 
4.325. 

The expected extinction in the regions relevant to the 
UZC survey is shown in Figure and the data set is 



shown in the same projection in Figure 20 for compari- 
son. In the northern sample, the extinction ranges from 
Am =0.012 to 1.7 with a median of 0.11. The cleanest 
spot in the southern sample has Am =0.061, and the ex- 
tinction gets as extreme as Am =63 in the Galactic plane 
(which we masked out). 

To get a first crude handle on the importance of ex- 
tinction, we applied our analysis separately to the rela- 
tively clean north and to the full, uncropped south (for 
this test, we used the regions delimited by solid lines in 
Figure |2(], not the dashed lines). Perhaps surprisingly, the 
South does not show a great power excess over the North 
even though the whopping extinction associated with the 
Galactic plane was included in this southern sample. This 
is presumably due to a combination of effects: The calcu- 
lations of Vogeley (1998) suggested that extinction (for the 
SDSS region in the north) would dominate over the cosmo- 
logical power spectrum only on scales k < 0.03//i^^Mpc, 
and the smaller volume of the UZC survey precludes us 
from effectively probing such large scales. Also, although 
extinction is severe near the Galactic plane, this is a rela- 
tively smooth feature and therefore does not greatly affect 
smaller scale fluctuations. 

Indeed, the power spectrum in the North is, if anything, 
somewhat greater than in the South. This agrees with the 
findings of Park et al. and da Costa et al. (1994), reflecting 
the fact that the nearby southern sky is more quiescent, 
lacking northern structures such as Virgo and the Great 
Wall. 

To be cautious, we nonetheless subjected the south- 
ern subsample to two additional cuts based on the dust 




Fig. 19. — The extinction Am predicted by Schlegel, J^'inkbeiner 
& Davis (1999) is shown in gnomonic equal-area projection with the 
North Celestial Pole in the center and right ascension o = at the 
top, increasing clockwise. The Galactic plane is seen to intersect 
parts of the southern survey region. Note that the South Galactic 
Hemisphere is at the top. 
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Fig. 20j^ The UZC galaxies are shown in the same coordinates 
as Figure hsi Note the visible deficit of galaxies in regions near the 
Galactic plane where Schlegel et al (1999) predict high extinction. 



map in Figure u3 before producing our main results, the 



power spectra shown in Section 4.1. Following Park et al. 
(1994), we moved the 4'* right ascension cut to 3''. We 
also excluded the region less than 13° from the Galactic 
plane. Since the north-south differences were relatively mi- 
nor without these cuts, we expect extinction to play only 
a subdominant role in this cropped data set. 

We also corrected the observed magnitudes with the 
Schlegel et al. map and considered creating a new magnitude- 
limited sample for analysis. Unfortunately, this would 
have reduced the galaxy count too drastically to be of 
much interest. Since the extinction gets as high as Am = 
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1.7 even in the north, the new magnitude cut would have 
to be shifted from 15.5 to 13.8 to ensure completeness, 
leaving only 8% of the galaxies from the original sample. 

4.5. Robustness to radial selection function errors 

To assess the extent to which radial selection function er- 
rors may be adding spurious power, we repeated our anal- 
ysis with a variety of different selection functions. We re- 
placed our Schecter parameters (—1.1,-19.3) by a grid in 
(a, M*)-space, always normalizing to match the observed 
number of galaxies. This had a substantial effect on the 
measured band powers, especially on the largest scales. 
We experimented with an iterative approach whereby the 
selection function was fine-tuned to minimize the large- 
scale power, but this unfortunately failed to converge: 
by grossly over-estimating the (small) number of distant 
galaxies, this procedure causes S = (n/n — 1) — s- —1, and 
since this function is spatially constant, it gives very little 
large-scale power. We therefore retained the radial se- 
lection function corresponding to the Schecter parameters 
measured from the data in the conventional way (de Lap- 
parent et al. 1989) in our quoted results above. However, 
in light of the sensitivity of the results to slight changes 
in the assumed luminosity function, the leftmost band- 
power measurement should be taken with a liberal dash 
of salt. In particular, the suspiciously high value seen at 
k ~ 0.03//i-iMpc for the 130/i-iMpc sample in Figure | 
may well be due to this effect. Selection function problems 
may be present even in the nominally volume-limited sur- 
veys, since Malmquist bias produces a selection function 
that is not quite constant near the far edge of the survey 
volume. 

5. CONCLUSION 

We have computed the redshift space power spectrum 
of the UZC galaxy redshift catalog. The results are sum- 
marized in Figures p[]lO| and are well fit by, e.g., a ACDM 
model with a moderate amount of bias (1.2 to 1.4). 

This paper is part of a larger effort to make galaxy red- 
shift survey analysis more comparable to the state of the 
art for CMB experiments, explicitly calculating the win- 
dow function for each band power measurement and in- 
cluding all finite- volume effects in the quoted error bars. 
However, it is merely the first step, and much work re- 
mains to be done: 

1. We measured only the power spectrum in redshift 
space. A more ambitious analysis should extract the 
additional information present in clustering anisotropics 
due to redshift-space distortions. 

2. We focused only on large scales, optimizing our 
pair weighting for the case of Gaussian fluctua- 
tions. For future work on smaller non-linear scales, 
k ^ 0.3//i~^Mpc, better results can be obtained 
by modeling the non-negligible correlations between 
different k-modes that have been observed in simu- 
lations (Mciksin & White 1999; Scoccimarro et al. 
1999; Hamilton 2000). 

3. We measured merely the galaxy power spectrum. A 
more ambitious analysis would attempt to model the 
biasing issues needed to translate this into a mea- 
surement of the underlying matter power spectrum. 



for instance by combining the stochastic biasing for- 
mahsm (Pen 1998; Dekel & Lahav 1999; Tegmark 
& Peebles 1998; Somerville et al. 2000) with red- 
shift space distortions as outlined in Tegmark (1998) 
and/or by studying higher-order moments. 

4. Given the sensitivity of the results to errors in the 
radial selection function, it would be worthwhile for 
a future analysis to include a likelihood analysis of 
its shape (Binggeh et al. 1988; Willmer 1997; Tresse 
1999). 

5. It is likely that model testing with future large sur- 
veys will require detailed Monte Carlo simulations 
to quantify and correct for subtle survey selection ef- 
fects {e.g., fiber collisions). Fortunately, the matrix- 
based analysis pipeline we have presented lends itself 
well to this: once the Fisher matrix has been com- 
puted, the time required to analyze another (Monte 
Carlo) survey scales merely as iVj, not as N'^. 

6. Our clustering analysis was limited to the power 
spectrum, i.e., to second moments of the density 
field. Higher-order moments of the density field 
are likely to contain a wealth of additional cosmo- 
logical information on small scales (see Frieman & 
Gaztaiiaga 1999 and references therein). 

These are major challenges, but the unprecedented quality 
of the impending data sets from ongoing redshift surveys 
such as 2dF and SDSS provide ample motivation to pursue 
them. 
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